#--------------------------------------------------------------------------------------------------------------------#
#----------------------------------------------- Portmann & Stojanovic ----------------------------------------------#
#--------------------------------------------------- November 2020 --------------------------------------------------#
#------- Are Immigrant-Origin Candidates Penalized in Virtue of Ingroup Favoritism or Outrgoup Hostility? -----------#
#--------------------------------------------------------------------------------------------------------------------#
#------------------------------- C Supplementary Tables and Figures for Main Analyses -------------------------------#
#--------------------------------------------------------------------------------------------------------------------#

rm(list=ls())
setwd(".../...")
load("d_nr15_cand")
load("d_nr15_cand_p")
load("d_nr15_cand_n")
load("d_nr15_cand_m")

# ----- Drop candidates on lists which are not classifiable on the left-right scale

d_nr15_cand <- d_nr15_cand[!d_nr15_cand$partypos == 4, ]


#-----------------------------------------------------------------------------------------------#
#---------------------------- Part 0: Models with fixed effects --------------------------------#
#-----------------------------------------------------------------------------------------------#

# ------------------------------------ Table 8 Online appendix: Estimating the relationship between candidate name (Swiss vs. non-Swiss) and negative preference votes

# ----- Estimate main models negative preference votes

m1 <- lm(npv_rel ~ nonswissname + bfscantno.x + list_id, data=d_nr15_cand_n)

m2 <- lm(npv_rel ~ nonswissname + rel_pos_c + bfsincumbency + precumulated + bfssex + age_c3 + partypos + bfscantno.x + list_id, data=d_nr15_cand_n)

m3 <- lm(npv_rel ~ nonswissname + partypos + rel_pos_c + bfsincumbency + precumulated + bfssex + age_c3 + partypos + partypos*nonswissname + bfscantno.x + list_id, data=d_nr15_cand_n)

# ----- Table

custom.name <- c("Intercept", "Non-Swiss name (0/1)", 
                 "Relative position on ballot", "Incumbent (0/1)",
                 "Pre-cumulated (0/1)", "Male (0/1)", "Age: 30-50 years", "Age: 50+ years", 
                 "Party position = center", "Party position = left",
                 "Non-Swiss name*Party position = center", 
                 "Non-Swiss name*Party position = left")

texreg(list(m1, m2, m3), booktabs=TRUE,
       caption="Estimating the relationship between a candidate's name (Swiss vs. non-Swiss) and negative preference votes",
       omit.coef = "bfscantno|id",
       custom.coef.names = custom.name,
       caption.above = TRUE, stars = c(0.001, 0.01, 0.05, 0.1))

# ------------------------------------ Table 9 Online appendix: Estimating the relationship between candidate name (Swiss vs. non-Swiss) and positive preference votes

# ----- Estimate main models positive preference votes

m1 <- lm(pos_mean ~ nonswissname + bfscantno.x + list_id, data=d_nr15_cand_m)

m2 <- lm(pos_mean ~ nonswissname + rel_pos_c + bfsincumbency + precumulated + bfssex + age_c3 + partypos + bfscantno.x + list_id, data=d_nr15_cand_m)

m3 <- lm(pos_mean ~ nonswissname + partypos + rel_pos_c + bfsincumbency + precumulated + bfssex + age_c3 + partypos + partypos*nonswissname + bfscantno.x + list_id, data=d_nr15_cand_m)

#num_listpos_c + 

# ----- Table

custom.name <- c("Intercept", "Non-Swiss name", 
                 "Relative list ranking", "Incumbent",
                 "Pre-cumulated", "Male", "Age: 30-50 years", "Age: 50+ years", 
                 "Party position = center", "Party position = left",
                 "Non-Swiss name*Party position = center", 
                 "Non-Swiss name*Party position = left")

texreg(list(m1, m2, m3), booktabs=TRUE,
       caption="Estimating the relationship between a candidate's name (Swiss vs. non-Swiss) and positive preference votes",
       omit.coef = "bfscantno|id",
       custom.coef.names = custom.name,
       caption.above = TRUE, stars = c(0.001, 0.01, 0.05, 0.1))

#-----------------------------------------------------------------------------------------------#
#----------------------------------- Part 1: Predicted values ----------------------------------#
#-----------------------------------------------------------------------------------------------#

# ----- Models

m1 <- lme(npv_rel ~ nonswissname + rel_pos_c + bfsincumbency + precumulated + bfssex + age_c3 + partypos + num_listpos_c + bfscantno.x, random = ~ 1 | list_id, data=d_nr15_cand_n)

m2 <- lme(pos_mean ~ nonswissname + rel_pos_c + bfsincumbency + precumulated + bfssex + age_c3 + partypos + num_listpos_c + bfscantno.x, random = ~ 1 | list_id, data=d_nr15_cand_m)

# ----- Predicted values negative preference votes

ef <- effect("nonswissname", m1)
summary(ef)
ef <- as.data.frame(ef)

# ------------------------------------ Table 10 Appendix: Predicted negative preference votes by candidate name (Swiss, non-Swiss)

ef$fit <- format(round(ef$fit, 3), nsmall = 3)
ef$se <- format(round(ef$se, 3), nsmall = 3)
ef$lower <- format(round(ef$lower, 3), nsmall = 3)
ef$upper <- format(round(ef$upper, 3), nsmall = 3)

print(xtable(as.matrix(ef), caption = "Predicted negative preference votes, by candidate name (Swiss, non-Swiss)"), include.rownames=FALSE)

# ------------------------------------ Figure 2 Appendix: Predicted values of negative and positive preference votes, by candidate name (Swiss, non-Swiss)

ef$fit <- as.numeric(as.character(ef$fit))
ef$lower <- as.numeric(as.character(ef$lower))
ef$upper <- as.numeric(as.character(ef$upper))

colors <- c("grey", "black")
pred_neg <- ggplot(ef, aes(nonswissname, fit, color=nonswissname)) + geom_point(color=colors, size=2) + 
geom_errorbar(aes(ymin=lower, ymax=upper), width=0.3, color=colors, size=1) + 
theme_bw(base_size=12) +
theme_bw() +
ylim(0.8, 1.06) +
labs(x="", y="Predicted negative preference votes"); pred_neg

#Predicted values positive preference votes

ef <- effect("nonswissname", m2)
summary(ef)
ef <- as.data.frame(ef)


ef$fit <- as.numeric(as.character(ef$fit))
ef$lower <- as.numeric(as.character(ef$lower))
ef$upper <- as.numeric(as.character(ef$upper))

# Figure
colors <- c("grey", "black")
pred_pos <- ggplot(ef, aes(nonswissname, fit, color=nonswissname)) + geom_point(color=colors, size=2) + 
geom_errorbar(aes(ymin=lower, ymax=upper), width=0.3, color=colors, size=1) + 
theme_bw(base_size=12) +
theme_bw() +
ylim(0.8, 1.06) +
labs(x="", y="Predicted positive preference votes"); pred_pos

pred <- ggarrange(pred_neg, pred_pos, 
labels = c("", ""),
ncol = 2, nrow = 1)
pred

png("pred.png", width = 8.5, height = 4.5, units = 'in', res = 300)
plot(pred)
dev.off()

# ------------------------------------------------------------------ Table 11 Appendix: Effect of candidate names on positive preference votes (models with Inverse gaussian distribution)

# ----- Models

m1 <- glmer(pos_mean ~ nonswissname + bfscantno.x +  (1|list_id), data=d_nr15_cand_p, family=inverse.gaussian(link="log"))

m2 <- glmer(pos_mean ~ nonswissname + rel_pos_c + bfsincumbency + precumulated + bfssex + age_c3 + partypos + num_listpos_c + bfscantno.x + (1|list_id), data=d_nr15_cand_p, family=inverse.gaussian(link="log"))

m3 <- glmer(pos_mean ~ nonswissname + rel_pos_c + bfsincumbency + precumulated + bfssex +  age_c3 + partypos + num_listpos_c  + bfscantno.x + partypos*nonswissname  + (1|list_id),  data=d_nr15_cand_p, family=inverse.gaussian(link="log"))

# ----- Table

custom.name <- c("Intercept", "Non-Swiss name", 
                 "Relative position on ballot", "Incumbent",
                 "Pre-cumulated", "Male", "Age: 30-50 years ", "Age: 50+ years", 
                 "Party position = centre", "Party position = left",
                 "Number of candidates on ballot",
                 "Non-Swiss name*Party position = centre", 
                 "Non-Swiss name*Party position = left")

texreg(list(m1, m2, m3), booktabs=TRUE,
       caption="Positive preference votes vs. non-Swiss name",
       omit.coef = "bfscantno",
       custom.coef.names = custom.name,
       reorder.coef = c(2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 1),
       caption.above = TRUE, stars = c(0.001, 0.01, 0.05, 0.1))

# ------------------------------------ Table 13 Appendix: Predicted positive preference votes by candidate name (Swiss, non-Swiss)

ef$fit <- format(round(ef$fit, 3), nsmall = 3)
ef$se <- format(round(ef$se, 3), nsmall = 3)
ef$lower <- format(round(ef$lower, 3), nsmall = 3)
ef$upper <- format(round(ef$upper, 3), nsmall = 3)

print(xtable(as.matrix(ef), caption = "Predicted positive preference votes, by candidate name (Swiss, non-Swiss)"), include.rownames=FALSE)

#-----------------------------------------------------------------------------------------------#
#---------------------------- Part 2: Distribution dependent variables -------------------------#
#-----------------------------------------------------------------------------------------------#

# ------------------------------------ Figure 3 Appendix: Distributions of dependent variables

#Distribution negative preference votes

norm <- rnorm(1000, mean = 1, sd = 0.25)
norm <- as.data.frame(norm)
dist_dv_neg <- ggplot(norm, aes(norm)) + geom_density(fill="grey40", alpha=0.5) + 
geom_density(data = d_nr15_cand_m, aes(npv_rel), fill="black", alpha=0.5) + 
theme_bw() + labs(x="Negative preference votes/Normal Distribution", y="Density"); dist_dv_neg

jpeg(file="dist_dv_neg.jpeg", 
width = 8, height = 5, units = "in", res = 300)
plot(dist_dv_neg)
dev.off()

#Distributuion positive preference votes (Inverse Gaussian Distribution shows good fit)

invgauss <- rinvgauss(1000, mean=3, shape=0.8)
invgauss <- as.data.frame(invgauss)
dist_dv_pos <- ggplot(invgauss, aes(invgauss)) + geom_density(fill="grey50", alpha=0.5) + 
geom_density(data = d_nr15_cand_m, aes(pos_mean), fill="black", alpha=0.5) + 
theme_bw() + 
labs(x="Positive preference votes/Inverse Gaussian Distribution", y="Density")

jpeg(file="dist_dv_pos.jpeg", 
width = 8, height = 5, units = "in", res = 300)
plot(dist_dv_pos)
dev.off()

#-----------------------------------------------------------------------------------------------#
#------------------ Part 3: Different specification of dependent variables ---------------------#
#-----------------------------------------------------------------------------------------------#

# --------- Models Cumulation and panachage as DV 

d_nr15_cand_c <- d_nr15_cand_m[!is.na(d_nr15_cand_m$cum_rel),]
d_nr15_cand_c <- d_nr15_cand_c[ ,!colSums(is.na(d_nr15_cand_c))]

d_nr15_cand_pan <- d_nr15_cand_m[!is.na(d_nr15_cand_m$pan_mean),]
d_nr15_cand_pan <- d_nr15_cand_pan[ ,!colSums(is.na(d_nr15_cand_pan))]

#Cumulation

con1 <- gamlss.control(n.cyc = 500)

m17 <- gamlss(cum_rel ~ nonswissname + bfscantno.x + re(random=~1|list_id), family=ZAIG, data=d_nr15_cand_c)

m18 <- gamlss(cum_rel ~ nonswissname + rel_pos + bfsincumbency + precumulated + bfssex + age_c3 + partypos + num_listpos_c + bfscantno.x + random(list_id, lambda = TRUE), family=ZAIG, data=d_nr15_cand_c, control=con1)

m19 <- gamlss(cum_rel ~ nonswissname + rel_pos + bfsincumbency + precumulated + bfssex + age_c3 + partypos + num_listpos_c + nonswissname:partypos + bfscantno.x + re(random=~1|list_id), family=ZAIG, data=d_nr15_cand_c)

#Panachage
m20 <- gamlss(pan_mean ~ nonswissname + re(random=~1|list_id), family=ZAIG, data=d_nr15_cand_pan)

m21 <- gamlss(pan_mean ~ nonswissname + rel_pos + bfsincumbency + precumulated + bfssex + age_c3 + partypos + num_listpos_c + bfscantno.x + re(random=~ 1| list_id), family=ZAIG, data=d_nr15_cand_pan, control=con1)

m22 <- gamlss(pan_mean ~ nonswissname + rel_pos + bfsincumbency + precumulated + bfssex + age_c3 + partypos + num_listpos_c + nonswissname:partypos + bfscantno.x + re(random=~1|list_id), family=ZAIG, data=d_nr15_cand_pan, control=con1)

# ------ Table 14 Appendix: Effect of candidate names (Swiss, non-Swiss) on votes from "cumulation" and "panachage"

custom.name.1 <- c("Intercept", "Non-Swiss name",
                   "Relative position on ballot", "Incumbent (0/1)",
                   "Pre-cumulated (0/1)", "Male (0/1)", "Age: 30-50 years ", 
                   "Age: 50+ years", "Party position = center", "Party position = left",
                   "Number of candidates on ballot","Sigma Intercept", "Nu Intercept")

texreg(list(m18, m21), booktabs=TRUE,
       caption="Effect of candidate name on panachage and cumulation",
       omit.coef = "bfscantno",
       custom.coef.names = custom.name.1,
       reorder.coef = c(2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 1, 12, 13),
       caption.above = TRUE, stars = c(0.01, 0.05, 0.1))


#-----------------------------------------------------------------------------------------------#
#--------------------- Part 4: different expalantory variables ---------------------------------#
#-----------------------------------------------------------------------------------------------#

#----------------------------- Analysis 1: Distinction name origins ----------------------------#

#Negative preference votes

m1 <- lme(npv_rel ~ origin_name_c3 + rel_pos_c + bfsincumbency + precumulated +  bfssex + age_c3 + partypos + num_listpos_c + bfscantno.x,  random = ~ 1 | list_id, data=d_nr15_cand_n)

m2 <- lme(npv_rel ~ origname_lang + rel_pos_c + bfsincumbency + precumulated + bfssex + age_c3 + partypos + num_listpos_c + bfscantno.x,  random = ~ 1 | list_id, data=d_nr15_cand_n)

# ------------------------------- Table 11 Appendix: Effect of candidate names (different outgroups) on negative preference votes

custom.name.1 <- c("Intercept", "Non-Swiss name, Western",
"Non-Swiss name, Non-Western", "Relative position on ballot", "Incumbent",
"Pre-cumulated", "Male", "Age: 30-50 years", "Age 50+ years", 
"Party position = center", "Party position = left",
"Number of candidates on ballot", "non-Swiss name from German-speaking country",
"non-Swiss name, from French-speaking country", "non-Swiss name, from Italian-speaking country", 
"non-Swiss name, from other country")

texreg(list(m1, m2), booktabs=TRUE,
caption="Effect of candidate names (different outgroups) on negative preference votes",
custom.coef.names = custom.name.1,
omit.coef = "cant",
caption.above = TRUE, stars = c(0.001, 0.01, 0.05, 0.1))

linearHypothesis(m1, "origin_name_c32=origin_name_c33") 
linearHypothesis(m2, "origname_lang1=origname_lang2")
linearHypothesis(m2, "origname_lang1=origname_lang3")
linearHypothesis(m2, "origname_lang1=origname_lang5")
linearHypothesis(m2, "origname_lang2=origname_lang3")
linearHypothesis(m2, "origname_lang2=origname_lang5")
linearHypothesis(m2, "origname_lang3=origname_lang5")


# ----------------------------- Table 15 Appendix: Effect of candidate names (different outgroups) on positive preference votes

m3 <- glmer(pos_mean ~ origin_name_c3 + rel_pos + bfsincumbency + precumulated +  bfssex + age_c3 + partypos + num_listpos_c + bfscantno.x + (1|list_id), data=d_nr15_cand_p, family=inverse.gaussian(link="log"))

m4 <- glmer(pos_mean ~ origname_lang + rel_pos + bfsincumbency + precumulated +  bfssex + age_c3 + partypos + num_listpos_c + (1|list_id), data=d_nr15_cand_p, family=inverse.gaussian(link="log"))

#Table

texreg(list(m3, m4), booktabs=TRUE,
caption="Effect of candidate names (different outgroups) on positive preference votes",
custom.coef.names = custom.name.1,
omit.coef = "cant",
caption.above = TRUE, stars = c(0.001, 0.01, 0.05, 0.1))

linearHypothesis(m3, "origin_name_c32=origin_name_c33") 
linearHypothesis(m4, "origname_lang1=origname_lang2")
linearHypothesis(m4, "origname_lang1=origname_lang3")
linearHypothesis(m4, "origname_lang1=origname_lang5")
linearHypothesis(m4, "origname_lang2=origname_lang3")
linearHypothesis(m4, "origname_lang2=origname_lang5")
linearHypothesis(m4, "origname_lang3=origname_lang5")




